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ABSTRACT 

We explore how radiative cooling, supernova feedback, cosmic rays and a new model of the energetic feed- 
back from active galactic nuclei (AGN) affect the thermal and kinetic Sunyaev-Zel'dovich (SZ) power spectra. 
To do this, we use a suite of hydrodynamical TreePM-SPH simulations of the cosmic web in large periodic 
boxes and tailored higher resolution simulations of individual galaxy clusters. Our AGN feedback simulations 
match the recent universal pressure profile and cluster mass scaling relations of the REXCESS X-ray cluster 
sample better than previous analytical or numerical approaches. For multipoles ( < 2000, our power spectra 
with and without enhanced feedback are similar, suggesting theoretical uncertainties over that range are rela- 
tively small, although current analytic and semi-analytic approaches overestimate this SZ power. We find the 
power at high 2000 - 10000 multipoles which ACT and SPT probe is sensitive to the feedback prescription, 
hence can constrain the theory of intracluster gas, in particular for the highly uncertain redshifts > 0.8. The 
apparent tension between cr 8 from primary cosmic microwave background power and from analytic SZ spectra 
inferred using ACT and SPT data is lessened with our AGN feedback spectra. 

Subject headings: Cosmic Microwave Background — Cosmology: Theory — Galaxies: Clusters: General — 
Large-Scale Structure of Universe — Black Hole Physics — Methods: Numerical 



1 . SZ POWER TEMPLATES AND THE OVERCOOLING PROBLEM 

When CMB photons are Compton-scattered by hot elec- 
trons, they gain energy, giving a spectral decrement in ther- 
modynamic tem£erature_bek2wv »_220 GHz, and an ex- 
cess above (ISunvaev & Zeldovichlll970l) . The high electron 
pressures in the intracluster medium (ICM) result in cluster 
gas dominating the effect. The integrated signal is propor- 
tional to the cluster thermal energy and the differential signal 
probes the pressure profile. The SZ sky is therefore an effec- 
tive tool for constraining the internal physics of clusters and 
cosmic parameters associated with the growth of structure, 
in particular the rms amplitude of the (lin ear) density power 
spectrum on cluster-m ass scales 0"s (e.g., iBirkinshawl 1 1 999t 
ICarlstrom et al J 120021) . Identifying clusters through blind SZ 
surveys and measuring the SZ power spectrum have been long 
term goals in CMB research, and are reaching f ruition through 
the South Pole Telescope, SPT dLueker et ail J 120091) and At- 
acama Cosmology Telescope, ACT dFowler et al] 120101) ex- 
periments. The ability to determine cosmological parame- 
ters from these SZ measurements is limited by the systematic 
uncertainty in theoretical modelling of the underlying cluster 
physics and hence of the SZ power spectrum. The power con- 
tribution due to the kinetic SZ (kSZ) effect that arises from 
ionized gas motions with respect to the CMB rest frame adds 
additional uncertainty. 

There are two main approaches to theoretical computations 
of the thermal SZ (tSZ) power spectrum: from hydrodynam- 
ical sim ulations of SZ sky maps or from semi-analytical es- 
timates (Honditan|200l |200l B0205). Large cosmolog- 
ical simulations providing a gastrophysical solution to the 
pressure distribution should include effects of non-virialized 
motions, accretion shocks, and deviations from spherical 
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symmetry. Averaging over many realizations of synthetic 
SZ sky projections yields the power spectrum and its vari- 
ance (e.g., B0205 ; ldaS7lvaetalll2000T: iSpringel etani200ll 
Seliaket al. 2001; Schaf er et alj|2006h . In conjunction with 
primary anisotropy signals and extragalactic source models, 
the SZ power spectrum has been used as a template with 
variable amplitude Asz for extracting cosmological parame- 
ters by the Cosmic Background Imager (CBI ) team (B0205; 
Sievers et all 120091) and the ACBAR team (Goldstei n et all 
2003t lReichardtet al.ll2009h . A sz was used to estimate a 
o"8,sz as a way to encode tension between the SZ- 

determined value and the (lower) erg obtained from the pri- 
mary anisotrop y signal. The CBI team also has included an 
analytic model (Ko matsu & Se liak 2002 , KS) which was als o 
the one adopted by the WMAP team dSpergel et al.ll2CK)7l) . 
The KS template yielded a lower value for <r 8 sz than that 
obtained with the simulation template, by ~ 10%. The KS 
model assumes a universal ICM pressure profile in hydro- 
static equilibrium with a polytropic (constant F) equation of 
state. The power spectrum is then obtained using an ana- 
lytic fit to 'halo model' abundances. So far the SPT and 
ACT have onl y used the KS template and a relate d semi- 
analytic oneJOsfr^eretaD 120051: iBode et al] 120091) . This 
model (Sehgal et al .1 l20l6T si0) allows map generation by 
painting dark matter halos in N-body simulations with gas. 
It expands on KS by calculating the gravitational potential 
from the DM particles, includes an effective infall pressure, 
adds simplified models for star formation, non-thermal pres- 
sure support and energy feedback which are calibrated to ob- 
servations. Using these templates, the SPT team derived a 
o~ r.sz lower than t he primary anisotropy <x 8 (e.g., WMAP7, 
Lars on et alJliOlOh . 

Current simulations with only radiative cooling and su- 
pernova feedback excessively over-cool cluster centers (e.g. 
iLewis et ail 120001) . leading to too many stars in the core, 
an unphysical rearrangement of the thermal and hydrody- 
namic structure, and problems when compared to obser- 
vations, in particular for the entropy and pressure pro- 
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files. The average ICM pressure profile found through 
X-ray observations of a sample of nearby galaxy clusters 
dArnaud et al.ll2009|) is inconsisten t with adaptive-mesh clus- 
ter sim ulations dNagai et al.ll2007l) . as well as the KS analytic 
model dKomatsu et al.ll2010l) . Pr e-heating (e.g . iBialek et al] 
2001 and AGN fee dback (e.g. ISiiacki et all 120071 12008b 



Puch wein et alJ 120081) help solve the over-cooling problem 



and improve agreement with observed c luster properties. 

Prev iously, an analytical model by Rovchowdhur v et all 
(2005) has explored the e ffects of effervescen t heating on 
the SZ power spectrum and iHolder et al.l (120071) use a semi- 
analytical model to calculate how an entropy floor affects the 
SZ power spectrum. There have been several simulations 
on galaxy and group scales that have studied how 'quasar' 
feedb a ck impacts the to t al SZ decrement dThacker et al j 
| 2006t Scannapieco etalj 120081: iBhattacharva et alJ l2008h 
Chatterjee et alJ I2008IK In this work we explore whether 
AGN feedback incorporated into hydrodynamical simulations 
of structure formation can suppress the over-cooling problem 
and resolve the current inconsistency between theoretical pre- 
dictions and observations of the SZ power spectrum and X-ray 
pressure profile. 

2. MODELED PHYSICS IN OUR SIMULATIONS 
2.1. Cosmological simulations 

We pursue two complementary approaches using smoothed 
particle hydrodynamic (SPH) simulations: large-scale peri- 
odic boxes provide us with the necessary statistics and vol- 
ume to measure the SZ power spectrum; individual cluster 
computations allow us to address over-cooling at higher reso- 
lution and compare our AGN feedback prescription with pre- 
vi ous models. We used a modified version of the GADGET- 
2 (Springel 2005) code. Our sequence of periodic boxes 
had sizes 100, 165, 330 hr x Mpc. The latter two used N DM = 
N g - ds = 256 3 and 512 3 , maintaining the same gas particle mass 
m gas = 3.2 x 10 9 h~ l M , DM particle mass m DM = 1.54 x 
lO lo /i _1 M and a minimum gravitational smoothing length 
s s = 20 h kpc; our SPH densities were computed with 32 
neighbours. For our standard calculations, we adopt a tilted 
ACDM cosmology, with total matter density (in units of the 
critical) Q m = Qdm + Qb = 0.25, baryon density £2b = 0.043, 
cosmological constant Qa = 0.75, Hubble parameter h = 0.72 
in units of 100 km s Mpc -1 , spectral index of the primor- 
dial p ower-spectrum n % _ .96 and cr 8 = 0.8. For the 'zoomed' 
cases (Katz & White 1993), we repeatedly simulated the clus- 
ter 'g676'(with the high resolution m gas = 1.7 X 10 8 /T 1 M , 
wdm = 1-13 x 10 9 h~ l M and e s = 5 h~_ kpc, using 48 neigh - 
bours to compute SPH densities, as in lPfrommer et al.ll 2007). 

We show results for three variants of gas heating: (1) 
the classic non-radiative 'adiabatic' case with only formation 
shock heating; (2) an extended radiative cooling case with star 
formation, supernova (SN) feedback and cosmic rays (CRs) 
from structure formation shocks; (3) AGN feedback in addi- 
tion to radiative cooling, star formation, and SN feedback. 
Radiative cooling and heating were computed assuming an 
optically thin gas of a pure hydrogen and helium primordial 
composition in a time-dependent, spatially uniform ultravi- 
olet background. Star formation and supernovae feedback 
were modelled using the hybrid multiphase mode l for the in- 
terstellar medium of Springel & Hernquist (2003a). The CR 
population is modelled as a relativistic population of protons 
described by an isotropic power-law distribution function in 
momentum space with a spectral index of a — 2.3, following 



lEnBlin et al.l ( 120071) . With those parameters, the CR pressure 
modifies the SZ effect at most at the percent level and causes 
a reduction of the res ulting integrated Compton-y parameter 
dPfrommer et alj|2007l) . 

2.2. AGN feedback model 

Current state-of-the-art cosmological simulations are still 
unable to span the large range of scales needed to resolve 
black hole accretion. Hence a compr omise treatment for 
AG N feedback is needed. F or example, ISiiacki et al.l d2007l) 
and B ooth & Schavd d2009l) adopted estimates of black hole 
accretion rates based on the Bondi-Hoyle-Lyttleton formula 
(Bondi & Ho vie 119441) . Here we introduce a sub-grid AGN 
feedback prescription for clusters that allows for lower res- 
olution still and hence can be applied to large-scale struc- 
ture simulations. We couple the black hole accretion rate to 
the global star formation rate (SF R) of the cluster, as sug- 
gested by Thompso n et al.l ([2005) using the following argu- 
ments. The typical black hole accretion rates and masses for 
the inner gravitationally stable AGN disks (of size < lpc) are 
~ lM /yr and ~ 1O 6 M . Since AGN lifetimes are much 
longer than 1 Myr, mass must be transferred from larger radii 
to the inner disk. However, at much larger radii this outer disk 
is gravitationally unstable and must be forming stars. Thus, 
in order to feed the AGN, stability arguments suggest that the 
rate of accretion must be greater than the SFR. For simplicity 
we assume that Mbh k M+. We inject energy into the ICM 
over a spherical region of size Ragn about the AGN, accord- 
ing to 



£inj(< #agn) = £ r AM< R AGN )c 2 At 

if M*(< R AGN ) > 5M /yr. 



(1) 



The duty cycle over which the AGN outputs energy is Af and 
e r is an 'efficiency parameter'. (As we describe below, the 
calculated efficiency for turning mass into energy is much 
smaller than e r .) We have explored a wide range of our two 
parameters, but the specific choices made for the figures are 
Af = 10 8 yr and e r = 2 x 10~ 4 . We require a minimum SFR of 
5 M /yr to activate AGN heating in the halo it is housed in. 

Given the output AGN energy, we must prescribe how it 
is to be distributed . Our procedure is motivated by the way 
ISiiacki & Springell ((2006) did AGN heating via bubbles. Us- 
ing an on-the-fly friends-of-friends (FOF) halo finding algo- 
rithm in GADGET-2, we determine the mass and center of 
mass of each halo with M ha io > 1.2 x 10 n h~ l M . We cal- 
culate its global SFR within the AGN sphere of influence of 
radius 



Ragn = max 



Mha 



do 
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xlOO/T 1 kpc, 
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where magn = £s and E(z) 2 - Q m (l + z) 3 + £2a- Within the 
halos we partition onto those gas particles inside of ^agn 
according to their mass. We have varied the prescription for 
Ragn and its floor magn (chosen here to be the gravitational 
softening e s ); the specific numbers given in eq. [2] (and for 
e r ) match previous successful model s that suppress the over- 
cooli ng b y means of AGN feedback (Sijacki et al. 200l|, see 
Sect. 13. It . Defining R A as the radius at which the mean inte- 
rior density equals A times the critical density p a (z) (e.g., for 
A = 200 or 500), then the ratio of Ragn to R200 is a constant 
~ 0.05. 
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Fig. 1 . — Shown are fi, (dashed lines) and f stM (solid lines) normalized to 
the universal value (fi, = 0.13) assumed in our simulations of our cluster 
g676 with M 500 = 6.8 X 10 13 h~ l M Q . The blue lines are for the simula- 
tion with radiative cooling and star formation while the red and orange lines 
are for our AGN feedb ack model (s r = 2 X 1CT 4 , M* S= 5M /yr) and that 
by ISiiack i et al. 12008), re spectively. The data points are observations by 
Gonzalez et al. 12007) and Afshordi et al. (2007). f stal (< #500) from X-ray 
measurements also agree s well, but the error s are large. Our sub-grid model 
matches the results from Siiacki et al. 12008) in this high resolution simula- 
tion well. 

Although we have referred to our feedback mechanism 
as being caused by AGN outflows, radiation pressure from 
stellar luminosity acting on dust grains will serve much the 
same purpose, and c ould also deliver high efficiencies (e.g. 
Thom pson et al.ll2~005l) . In the code, we have so far added 
as a pure heating component, but it should allow for a me- 
chanical, momentum-driven wind component as well, which 
would not be as prone to catastrophic cooling and likely de- 
crease the s r needed for useful star formation suppression. 

The relevant energy budget is not in fact defined by e r , but 
rather by a redshift-dependent effective feedback efficiency 
£ e ff = £,'£'inj,i7[M*(< r) c 2 ], where we sum over every energy 
injection event (labeled by i) and we calculate the stellar mass 
M+(< r) within a given radius. In all cases, e e ff <K £ r , be- 
cause: (i) heating suppresses the stellar mass AM* created 
over At, making it quite a bit less than the stellar mass M*Af 
that would have formed without any feedback; and (ii) E m j is a 
stochastic variable, which we find to be zero about half of the 
time because the required SF threshold is not achieved. With 
our fixed e, - Racn prescription, our canonical g676 example 
has e e n" ~ 5 x 10~ 6 for the entire simulation; if all energy had 
been released within the final Ragn, £eff would be 8 X 10" 5 , 
but feedback, especially at early times, is much more widely 
distributed. Of a total E mi = 9 x 10 61 ergs for g676 we find 
58% is delivered in the cluster formation phases at z > 2, an- 
other 23% is delivered in the redshift range 1 < z < 2 that 
can be probed with ACT and SPT resolution, and only 19% 
comes from the longer period below redshift 1 . Feedback pre- 
scriptions with smaller which still give the desired star 
formation suppression need further exploration. 

3. PRESSURE PROFILES 

3.1. Testing AGN feedback as resolution varies 

AGN feedback self-regulates the star formation and ener- 
getics of a cluster. In Fig. Q] we compare the fraction of 
baryons (/b) and stars (/ star ) as functions of cluster radius for 



the high-resolution 'g676' simulations. Our radiative simula- 
tion produces 1.5-2 times more stars than those with AGN 
feedbac k. Our sub-grid AG N model nicely reproduces the re- 
sults in Sii acki etaLI d2008h . It should also produce reliable 
results in the cosmological box simulations in which over- 
cooling is less severe because of the lower resolution. There 
is significant sensitivity to the value chosen for the feedback 
parameter e r : doubling it lowers fi, by a factor of 1.5, halv- 
ing it increases / stal - by 1.4. The 100 h Mpc simulations 
were used to study the resolution dependence of our feed- 

1 13 

back model by varying A^ gas in steps from 64 to 256, with e s 
and hence magn (eq. [2 decreased accordingly. As magn de- 
creased, / stal within /?5oo increased almost linearly for radia- 
tive cooling, whereas for AGN feedback the increases were 
less. This can be traced to the hierarchical growth of structure 
since in low-resolution simulations: the small star forming 
systems are under-resolved; this decreases the SFR that me- 
diates our AGN feedback; and this lowers the overall number 
of stars produced in the simulati ons. This behaviou r is seen 
in other AGN feedback models (Siiacki et al. 2007) and has 
be en extensively studied in non -AGN feedback simulations 
bv lSpringel & Hernquistl(l2003bl) . 

3.2. Stacked pressure profiles 

For every halo identified by our FOF algorithm, we calcu- 
late the center of mass, R A , the mass Ma within R A and com- 
pute the spherically-averaged pressure profil e normalize d to 
P A = GM A Ap CI (z)f b /R A , with f b = Q b /Q m (lVoitl l2005h and 
radii scaled by R A . We then form a weighted average of these 
profiles for the entire sample of clusters at a given redshift. 
For Fig. |2] we have weighted by the integrated y-parameter, 

cr C Ra 

Y A = —^- P e (r)4nr 2 dr cc £ th (< fl A ) , (3) 
(m e c z ) Jo 

where <x x is the Thompson cross-section, m e is the electron 
mass and P e is electron pressure. For a fully ionized medium 
the thermal pressure P = P e (5X H + 3)/2(X H + 1) = 1.932P e , 
where Xu = 0.76 is the primordial hydrogen mass fraction. 
Splitting the clusters into a number of mass bins gives similar 
results to this monolithic Y A weight, as does weighting by Y^. 
We have found that a simple parametrized model 

P/P 500 = A [1 + (x/x c rr /a , x = r/R 5Q0 , (4) 

with core-scale x c , amplitude A, and two power law indices, 
a and y, fits better than with a fixed a. Sample values for our 
AGN feedback are A = 82, x c = 0.37, a = 0.84 and y = 4.6 
at z — 0; generally the parameters depend upon cluster mass 
and redshift. At z > 1, a more complex parameterization is 
needed. 

In Fig. |2] we show average pressure profiles multiplied by 
x 3 to make them oc dE^/dlnr, the thermal energy per loga- 
rithmic interval in radius, and hence to dY A /dlnr. All pro- 
files of dE t h/dlnr from simulations and observations peak 
at or before P200, but an integration to at least 4/?2oo is re- 
quired for the total thermal energy to converge. By contrast, 
the KS profile does not drop over this range due to the con- 
stancy of F and does not include the outer cluster phenom- 
ena of asphericity, accretion shocks, etc. Throughout this pa- 
per, we have computed the K S model with an up dated con- 
centration parameter given by iDuffv et al] (120081) . We also 
show a scaled average S10 pressure profile for clusters with 
1O 14 M < M 5Q o < 5xl0 14 and redshift < 0.2. The S 10 profile 
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Fig. 2. — Top: Comparison of fits to normalized average pressure profiles 
from analytic calculations, simulations and observations, scaled by (r/Rsoo) 3 - 
For a cluster of M500 = 2 X 10 14 /j~'M o , we show fits to our SPH sim- 
ulations (red), and compare them with the analytic KS profile (green), the 
semi-analytic S10 average profile (light green), and a fit to AMR simulations 
(updated profile by Nagai et al. 2007, private communication; orange). Our 
feedback model matches a fit to X-ray observations ( Arnaud et al. 2009, grey 
bands) within R500 well; only the dark grey part is actually a fit to the data, 
with the light grey their extrapolation using older theory results unrelated to 
the data. We illustrate the 1 and 2 cr contributions to centered on the me- 
dian for the feedback simulation by horizontal purple and pink error bars. 2nd 
panel: We compare fits to our AGN model at redshift z = (red solid) to all 
our three models at redshift z = 1 (blue). Shown are the lcr error bars of the 
cluster-by-cluster variance of the weighted averages in our AGN models us- 
ing corresponding lighter colors. 3rd panel: We show the effective adiabatic 
index T for our simulations, comparing it with KS (dash-dotted) and with a 
constant 1.2 (light green). Bottom: The distribution of kinetic-to-thermal en- 
ergy in percentile decades is indicated by the dots for the feedback case, with 
the median shown for all three models; thus, there are significant additions to 
pressure support even in the cores of simulated clusters, and even more so in 
the SZ-significant outer parts. 

has been weighted by and agrees well within R500 and with 
a slight excess pressure beyond R500. 

Fig.[2]shows our fee dback model t r aces t he observed "uni- 
versal" X-ray profile of lArnaud et al.l (|2009) shown as a dark- 
grey band rather well within R500. This fit came out natu- 
rally, with no further tuning o f our feedback par ameters be- 
yond trying to agree with the Siiac ki et al.l ((2008) simulation. 



Our models without AGN feedback have larger pressures in- 
side /?3oo. For the light grey band beyond the universal 
X-ray profile did not use observations, but was fit to an av- 
erage profile of earlier simulations so the deviation > R500 
does not represent a conflict of our profiles with the data, 
rather with the earlier simulations. The band shown for the 
X-ray profile gives a c rude correction for th e bias in M500 and 
7?5oo resulting from the Arnau det al.l (120091) assumption of hy- 
drostatic equilibrium. This yields mass values which are on 
average 25% too low dNagai et al.ll2007l) . so the band repre- 
sents a 0-25% uncertainty in M50Q. This change only affects 

R500 oc M™ and Psoo^oq rc ^500 but does not a ff ect the 
shape of the profile. (However, as the bottom panel shows, 
such a correction fro m turbulence and un-virialized bulk mo- 
tions dKravtsov et al J 120061) will depend upon radius and se- 
lection function of the X-ray clusters used to make the fit.) 

Another important issue is the relation between the and 
cluster mass. We fit our results for this to the scaling relation 

T 50 o = 10 B (M500/3 x 10 14 h- x M G ) A h'l 12 Mpc 2 , (5) 

where hjo = 0.7 x 100 km s Mpc -1 . Parameters from our 
simulation are B = (-4.47 ± 0.08, -4.46 + 0.20, -4.5 + 0.1) 
and A = (1.66 + 0.12, 1.71 + 0.25, 1.75 + 0.06) for the se- 
quence (1) shock heating, (2) radiative cooling and (3) AGN 
feedback. These values are simi lar to the B = -4.739 + 0.003 
andA= 1 .790 ±0.015 found bv lArnaud et all (120091). as well 
as the B = -4.713 + 0.004 and A = 1 668 + 0.009 found by 
S10. We note that Arna ud et al.l ((2009) actually used a mass 
proxy in place of M500, so their errors are not representative of 
the true observational sca tter in the Y—M scal ing relation. The 
AGN feedback model of Siia cki et all d2007l) was also able to 
reconcile th e cluster X-ray lumino sity and temperature scal- 
ing relation dPuchwein et al. 2008). 

We find a large variation in the outer pressure profiles be- 
yond R v i r , especially at redshift z ~ 1 as is shown in the sec- 
ond panel of Fig. [2] These regions may have sub-halos, and 
external but nearby groups on filaments, most of which will 
eventually be drawn into the clusters. In spite of the large 
variance of the scaled profiles, the fit to the profiles at z — 
follows the average. At larger redshift, however, our fitting 
formula will require more degrees of freedom than in eq.|5]to 
reflect the range of behaviour of the highly dynamical outer 
regions. 

4. SZ POWER SPECTRA FROM HYDRODYNAMICAL SIMULATIONS 

4. 1 . Stacked SZ power spectra of translated- rotated 
cosmological boxes 

We randomly rotate a nd translate our simulation snap- 
shots at different redshifts (Ida Silva et al.ll2000l:ISpringel et ail 
2001, B0205). To obtain thermal Compton-y maps, we per- 
form a line-of-sight integration of the electron pressure within 
a given solid angle, i.e. y — <x x J n e kT e /(m e c 2 )dl, where 
k is the Boltzmann constant, n e and T e are the number den- 
sity and temperature, respectively. We construct 1 .6°x 1 .6°and 
3.2°x3.2°maps for the 256 3 and 512 3 simulations, respec- 
tively. Using this method there are large sample variances 
(White et al. 2002) associated with nearby cluster contamina- 
tion. We have quantified their influence on the power spec- 
trum for each of our three physics models by averaging over 
twelve translate-rotate viewing angles each projected from 
our ten 256 3 full hydrodynamical simulations for each of the 
33 redshift outputs back to a redshift z = 5; the power spectra 
of which are then added up to yield the total spectrum. This 
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Fig. 3. — Predictions for the tSZ power spectrum a t 30 GHz fr om our simula- 
tions (red and purple symbols), simulations by Springel et al. ( 2001) (orange 
triangles), simulations by Bond et al. 1 2005) (orange pluses), semi-analytical 
simulations by S10 (dark green) and analytical calculations by KS (light 
green). The 256 3 power spectra (red symbols) are averages over 12 translate- 
rotate tSZ maps and 10 separate hydrodynamical simulations for each of the 
33 redshift bins, the power spectra of which are then added up to yield the 
total spectrum; the error bars show the variance among the power in all maps. 
The full-width half-max values appropriate for Planck, ACT and SPT show 
which part of the templates these experiments are sensitive to. At low-f, the 
discrepant higher power in the semi-analytical calculations can be traced to 
the enhanced pressure structures assumed beyond S200 over what we find. 

method of computing the power spectrum has the advantage 
of taking care of the artificial correlations that occur because 
any individual simulation follows the time evolution of the 
same structure. For the shock heating case, we did ten more 
hydrodynamical simulations to show that our averaged tem- 
plate had converged (within ~10%), but note that using only 
a few boxes can be misleading in terms of rare events. 

The computationally more expensive 512 3 SZ spectra have 
the equivalent of 8 256 3 plus wider coverage, so the 512 3 
shock heating result shown gives a reasonable indication of 
what to expect. The other 2 physics single-box cases at 512 3 
are similar to the 256 3 ensemble means. The analytical ap- 
proach has the great advantage of including an accurate mean 
cluster density to high halo masses, but to be usable for SZ 
power estimation, scaled pressure profiles must also be accu- 
rate, a subject we turn to in future work. For now, we note that 
using such profiles from our simulations gives good agree- 
ment with the average SZ power shown at the low t where 
sample variance will be largest. In Fig. [3] our simulation tem- 
plates and the KS template shown have excluded structures 
below z = 0.07 to decrease the large sample variance associ- 
ated with whether a large-ish cluster enters the field-of-view. 
Such entities would typically be removed from CMB fields 
and considered separately. 

The mean Compton >'-parameter found in our AGN feed- 
back simulations is one order of magnitude below the COBE 
FIR AS upper limit of 15 x 10" 6 dFixsen et al.ll 19961) . 

We compare the theoretical predictions for the tSZ power 
spectrum in Fig. [3] Our 512 3 and 256 3 shock heating 
simulations are in agreement with p revious SPH simulation 
power spectra ( Springel et al. 2001, B0205) scaled by Q oc 
(Q.\,h) 2 Q. m (Tg, with the factors determined from our simula- 
tions of differing cosmologies. The B0205 SZ power shown 
had a cut at z = 0-2, appropriate for CBI fields; using the same 
cut on a shock heating simulation with the same cosmology 



that we have done, we get superb agreement. 

The KS and S10 semi-analytic SZ power spectra templates 
differ substantially from our templates, in particular with 
higher power at low I: as shown in Fig. |2j the KS pressure 
profile beyond overestimates the pressure relative to both 
simulations and observations, leading to the modifie d shape 
and la rger this behaviour is also shown in Komat su et alj 
(2010). The spectrum from S10 is very similar to KS possi- 
bly because both assume hydrostatic equilibrium, and a poly- 
tropic equation of state with a fixed adiabatic index, T ~ 
1.1 - 1.2. Inside /?200, these assumptions are approximately 
correct, but they start to fail beyond /?200- A demonstration 
of this is the rising of F and of the ratio of kinetic-to-thermal 
energy K/U shown for our simulations in the bottom panels 
of Fig. [3] The present day (a — 1) internal kinetic energy 
of a cluster is given by K = E;m gaSj j \v-, - v + Ho(x-, - x)\ 2 /2, 
where Hq is the present day Hubble constant, Vj and ;t; are 
the peculiar velocity and comoving position for particle i, and 
v and x are the gas-particle-averaged bulk flow and center 
of mass of the cluster. The additional thermal pressure sup- 
port we find at large radii from AGN feedback results in the 
slightly slower rate of K/U growth shown. In all cases the 
large kinetic contribution shown should be properly treated in 
future semi-analytic models. 

Varying the physics over the three cases for energy injec- 
tion in our simulations leads to relatively minor differences 
in Fig. [3] among the power spectra for ( < 2000. This agree- 
ment is due in part to hydrostatic readjustment of the structure 
so the virial relation holds, which relates the thermal content, 
hence Fa, to the gravitational energy, which is dominated by 
the dark matter. Our AGN feedback parameters do not lead to 
dramatic gas expulsions to upset this simple reasoning. Our 
radiative cooling template has less power at all scales com- 
pared to the shock heating template since baryons are con- 
verted into stars predominantly at the cluster centers and the 
ICM adjusts adiabatically to this change. Thus, at low I where 
clusters are unresolved, shock heating and radiative simula- 
tions give upper and lower limits, bracketing the AGN feed- 
back case. AGN feedback suppresses the core value of the 
pressure compared to the radiative simulation resulting in less 
power at i > 2000, a trend that is more pronounced at z > 1 
(as shown in Fig. [3}. Thus, at these angular scales, the power 
spectrum probes the shape of the average pressure profile. It 
depends sensitively on the physics of star and galaxy forma- 
tion e.g..lScannapiec o~et aljj200 8). Over the Grange covered 
by Planck, these effects are sub-dominant, and serve to high- 
light the importance of the high-resolution reached by ACT 
and SPT 

4.2. Current constraints on SZ template amplitudes and erg sz 

Instead of varying all cosmological parameters on which 
the thermal and kinetic SZ power spectra, Q jt sz and Q ; ksz, 
depend, we freeze the shapes by adopting the parameters for 
our fiducial cr 8 = 0.8 (and Q^h = 0.03096) model evaluated 
at 150 GHz, and content ourselves with determining template 
amplitudes, A t sz and Aksz, and a total SZ amplitude Asz: 

AszQ.sz = /(v)A t szCV jt sz + ^kszQ.ksz ■ (6) 

The spectral function for the tSZ (ISunyaev & Zeldovichl 
119701) . /(v), vanishes at the SZ null at ~ 220 GHz and we 
normalize it to unity at v = 150 GHz, so it rises to 4 at 30 
GHz. Therefore if we find values of Asz below unity then 
either <x 8 is smaller than the fiducial cosmological value as 
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Fig. 4. — Our 150 GHz tSZ adiabatic and feedback (Asz = 1) power spectra 
computed with <x 8 = 0.8 (long dashed lines) are contrasted with the dark grey 
band indicating the lrr range in multiplicative amplitude, Asz = 0.75 ± 0.36, 
allowed by the SPT__FG power spectrum for the feedback template shape. 
The light grey band is the 2<r upper limit region. The Asz = 1 S10 tSZ power 
spectrum (dashed line) and the KS tSZ spectrum (dash dotted line) are shown 
for contrast; their allowed lcr band is determined by multiplying these by 
their Asz values given in Table 1, but cover a similar swath to the grey bands. 
We also show the averaged kSZ power spectra computed for our simulations 
by dotted lines. The kSZ spectra were calculated in the same was as the tSZ 
spectra were, and have similar shapes. However, kSZ is underestimated at 
low C because of missing bulk velocities in the simulations. There should 
be an additional (rather uncertain) kSZ template from inhomogeneous re- 
ionization as well. To show the tension with the CMB data, we plot the tSZ 
+ 0.46 kSZ power (solid lines) since this can be directly compared with the 
SPTdsfg grey bands. 

derived from the primary CMB anisotropics, or else the theo- 
retical templates overestimate the SZ signal. 

To determine the probability distributions of these am- 
plitudes and other cosmological parameters from cur- 
rent CMB data we adopt Markov Chain Monte Carlo 
(MCMC) techniques u sing a modified version of CosmoMC 
.Lewis & Bridld 120021) . We i nclude WMAP7 .Larson et al.l 
120101) and, separat ely, ACT .Fowler et al.l 120101) and SPT 
dLueker etafl l2009). In all cases, we assume spatial flatness 
and fit for 6 basic cosmological parameters (Clbh 2 , ^dm/. 2 * 
n s , the primordial scalar power spectrum amplitude A s , the 
Compton depth to re-ionization r, and the angular parame- 
ter characterizing the sound crossing distance at recombina- 
tion 6). We also allow for a flat white noise template C^ src 
with amplitude A sre , such as would arise from populations 
of unresolved point sources. We marginalize over A src , al- 
lowing for arbitrary (positive) values. Generally there will 
also be a spatial clustering component for such sources, and 
these will have templates that are partially degenerate in shape 
with that for tSZ, but because of the large uncertainties we 
ignore such contributions here. Reducing the SZ and unre- 
solved source problems to determinations of overall ampli- 
tudes multiplying shapes has a long history, e.g., the CBI se- 
quence of papers, and was adopted as well by the ACT and 
SPT teams. Our results differ slightly from those reported by 
the ACT team because they use WMAP5+ACT and a com- 
bined tSZ+kSZ SlO-template, and by the SPT team who use 
WMAP5+QUaD+ACBAR+SPT and add constraints on the 
white noise source amplitude beyond the non-negativity we 
impose. 

We first consider a simplified case with A_sz constrained to 
be zero and all other cosmic parameters and the source ampli- 



tude marginalized, yielding a probability distribution for Asz- 
The means and standard deviations from our MCMC runs are 
given in the upper rows of Table [T] in columns 2, 4 and 6 for 
a number of data combinations and for our 3 physics simu- 
lation cases, contrasting with KS and S10. The ACT data is 
for 148 GHz. There are two SPT cases given. The first uses 
just the 153 G Hz spectrum so it ca n be directly compared to 
ACT. For SPT, iLueker et al.l ((2009) also report a power spec- 
trum derived from subtracting a fraction x of their 220 GHz 
data from the 153 GHz data to minimize the contribution from 
dusty star-forming galaxies (DSFG); since 220 GHz is the SZ 
null, this does not modify the tSZ contribution, but would di- 
minish the frequency-flat kSZ. However, a normalization fac- 
tor is chosen to preserve power for primary CMB signals that 
are flat in frequency like kSZ. This has the effect of bo osting 
the tSZ power by a factor of ( 1 - x)~ 2 . ILueker et ail (120091) find 
that x = 0.325 minimizes the contribution from the DSFGs so 
the DSFG-subtracted spectrum suppresses the kSZ by a fac- 
tor of 0.46 relative to the tSZ. A ~25% uncertainty remains in 
x which should be taken into account statistically, but is not 
here. The correct approach would be to simultaneously treat 
the 153 GHz and 220 GHz cases, with full modelling of the 
different classes of point sources, including their clustering, 
and to take into account the non-Gaussian nature of the SZ 
and source signals which impact sample variance. 

The ACT data is only giving upper limits with their cur- 
rent published data, whereas SPT has detections at 153 GHz 
with Asz compatible with unity. For the SPT 153 GHz-only 
spectrum, we find S10 gives Asz = 1.39 + 0.34 while the 
feedback template gives Asz = 1-76 + 0.43, and the compara- 
ble 95% upper limits from ACT are 1.95 and 2.93. However, 
although the white noise shape has been vetoed by marginal- 
ization, there could be a residual clustered source contribution 
from dusty galaxies pushing the derived Asz high- To the ex- 
tent that SPT_)_ pc vetoes this DSFG clustering as well as their 
Poisson contribution, that Asz would be a better indicator. It 
shifts from 0.43 + 0.21 for KS and 0.50 + 0.25 for S10 up to 
0.75 + 0.36 for the feedback template, an increase of 50%. 
The large difference between the 150 and source-subtracted 
templates, even after marginalizing over a Poisson term, may 
suggest the power in the correlated source component may be 
similar to the SZ power, emphasizing the work necessary to 
do a correct treatment. 

Any non-zero kSZ contribution will take some of the am- 
plitude from Asz, leaving even smaller A t sz values; columns 
3, 5 and 7 of the table give estimates of this diminution. The 
kSZ power spectra that we have computed are broadly similar 
to the tSZ power shape, with however sufficiently significant 
differences to allow shape discrimination in addition to the 
frequency separability, as Fig. |4] shows. At 150 GHz and an 
I = 3000 pivot, we find the kSZ power is ~ 29%, ~ 29% and 
~ 27% of the tSZ power for the shock heating, radiative cool- 
ing and feedback simulations, respectively. We normalize the 
kSZ to the tSZ at this pivot of 3000 since it has most of the 
constraining power in the CosmoMC chains for the ACT and 
SPT measurements and results in the smallest error bars: on 
larger scales, the errors are increased by the contribution from 
primary anisotropics while smaller scales are dominated by 
the instrumental and galaxy-source shot noise. 

We used exactly the same procedure to obtain the kSZ 
spectrum as we used for the tSZ spectrum. The temperature 
decrement due to the kSZ effect is AT/T - crj J n e vjcdl, 
where v r is the radial peculiar velocity of the gas relative to 
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TABLE 1 

Constraints on Asz and 0"8.sz 



tSZ template 


ACT 148 GHz 


SPT 153 GHz 


SPTdsfg 




Asz 


A,sz 


Asz 


Asz 


Asz 


A t sz 


KS 


< 1.55 


< 1.26 


1.01 ±0.25 


0.72 ± 0.25 


0.43 ± 0.21 


0.30 ±0.21 


SIO 


< 1.95 


< 1.67 


1.39 ±0.34 


1.11 ±0.34 


0.50 ± 0.25 


0.38 ± 0.25 


Shock heating 


< 2.13 


< 1.84 


1.13 ±0.28 


0.84 ± 0.28 


0.44 ± 0.22 


0.31 ±0.22 


Radiative cooling 


<2.75 


< 2.46 


1.50 ±0.37 


1.21 ±0.37 


0.59 ± 0.29 


0.45 ± 0.29 


Feedback 


<2.93 


<2.66 


1.76 ±0.43 


1.49 ±0.43 


0.75 ± 0.36 


0.63 ± 0.36 




0"8,SZ 


0"8,tSZ 


0"8,SZ 


0"8,ISZ 


0"8,SZ 


0"8,ISZ 


KS 
SIO 

Shock heating 
Radiative cooling 
Feedback 


< 0.864 

< 0.891 

< 0.900 

< 0.935 

< 0.944 


< 0.845 

< 0.874 

< 0.883 

< 0.922 

< 0.932 


n 709+0.029 

u - ' y -0.029 
n Q9Q+0.031 
U.828_ 0030 

0.804« 

0.837*"? 
0.856^3 


fl 7S7+ 039 
u -' J '-0.037 

0.800!;| 
0.768!;g 
809+° 039 

O.OO7_ 0037 

n s^s+0.038 

U - SJ5 -0.037 


690+ 057 

-0.055 

705 +0 06 ° 

u -' -0.058 
06gl +0.059 
v m -0.057 
n 79 1 +0.060 
u -' -0.058 

7 46+ 061 

-0.055 


622+ 0105 

u - u -0.051 

636+ 0109 

-0.054 

607+° 119 
u,ou -0.057 

660+° 102 

-0.053 

703+° 091 

-0.055 



The mean and standard deviation of the thermal SZ power spectrum template amplitude A t sz and the total SZ, including our computed kSZ 
contribution. The numbers assume the kSZ template is perfectly degenerate in shape with the tSZ one. Asz = A t sz + A^sz at 150 GHz, with the 
relative enhancement in our simulations given by A^sz/Atsz = 0.29, 0.29, 0.27 for the shock heating, radiative cooling and feedback simulations, 
respectively. We have used the ACT team's 148 GHz power spectrum, the SPT team's 153 GHz spectrum and the SPT DSFG-subtracted 
(SPTosfo) spectrum, along with WMAP7. The amplitude of the SZ power is normalized to our fiducial erg = 0.8 cosmology. A rough guide to 
the o"8 tension is obtained in the lower rows, using crgsz 00 A ' (fl^A) -2 ' 7 , with exponents determined by B0205 and KS. Since kSZ varies more 
slowly with <xg than tSZ, the numbers are just indicative. 



the observer. We constructed 12 translate-rotate kSZ maps 
for each of our 10 separate hydrodynamical simulations and 
for each of the 41 redshift bins back to z = 10 (rather than 
Z = 5 for tSZ), computing the average and variance of all 
of these. Since we use simulations with side length L = 
165 /r'Mpc for our 256 3 cases, with fundamental wavenum- 
ber (26 /i^'Mpc) -1 , our spectra are missing a bit of power on 
the largest scales (affecting low-^) since we do not sample 
well the long-wavelength tail of the velocity power spectrum 
in spite of the number of runs done. 

We have included the kSZ template by ignoring the rela- 
tively small shape difference about the pivot point of the kSZ 
and tSZ power spectra; i.e., we assume the perfect degeneracy 
Cf.ksz = Q.tsz, as the SPT team did. Thus we only need the 
ratios AkszAAsz given above for the 150 GHz cases and the 
further x factors for the mixed frequency DSFG case. For the 
ratios we use our translate-rotate values of 0.29, 0.29 and 0.27 
from our simulations, 0.276 for S10, and used a rough esti- 
mate of 0.25 for KS. Apart from ignoring the shape difference, 
we have also ignored kSZ from patchy re-ionization at high 
redshift, although it can have a competitive amplitude to the 
late time fully ionized gas motion s with respect to the C MB 
rest frame that we are modelling dlliev et alj|2008l l2007h . In 
presenting the results from our analyses of the MC Markov 
chains, we just subtract A^sz from Asz- The Table[T]A t sz that 
we derive from these assumptions are all on the low side of 
unity for DSFG, with KS and S10 being more than 2.5cr low, 
whereas the feedback template is only about lcr low (and l<x 
high for 153 GHz alone). We leave it to future work to include 
a more complete implementation of the kSZ spectra. 

The means and errors on Asz provide the cleanest way of 
presenting the tension, or lack thereof, of these SZ models 
with the primary CMB data which indicates <xg ss 0.8. How- 
ever, it has been conventional to translate these numbers into 
a o"8,sz using the way Asz scales with cosmic parameters, 
roughly as Asz K 0"g (Q.\,h) , as given by B0205 and KS. The 
lower rows in Table 1 show crgsz using this scaling. Although 
the scaling applies to the tSZ component only, with the kSZ 
power being less sensitive to o~s, we also quote results for the 



kSZ-corrected cases. Ideally one should use the data to deter- 
mine the cosmic parameters which uniquely and fully deter- 
mine the primary spectrum, the A t sz and Aksz, and the tSZ and 
kSZ shape modifications as the parameters vary. This slaved 
treatment enforcing cr 8 S z = cr 8 has erg's value being driven 
by WMAP7 and other primary CMB data rather than by the 
SZ information. 

5. CONCLUSIONS AND OUTLOOK 

Without hydrodynamical simulations in a cosmological 
framework similar to the ones presented in this paper it is hard 
to come up with a consistent model of the gas distribution in 
clusters and the infall regions which both contribute signifi- 
cantly to the SZ power spectrum. In this paper, we identify 
three main points that a future semi-analytic model of such a 
pressure distribution has to provide. 

(1) In order to arrive at a consistent gas distribution that 
matches not only the integrated stellar mass fraction but also 
the X-ray derived pressure profiles within R500, we need self- 
regulating AGN-type feedback. We emphasize that we tuned 
our parameters to match a previous single-cluster model that 
successfull y suppressed the o ver-cooling by means of AGN 
feedback (Siiacki et al. 2008). The excellent agreement with 
current data was a pleasant byproduct: our simulated pres- 
sure profiles agree with recently obtained observational ones 
that have been constructed from X-ray data; the scaling rela- 
tions between the cl uster mass and X-ray based Compton-T 
dArnaud et all 12009) also agree; as do the integrated stellar 
and g as mass fractions ([Gonzalez et al. 2007; Afshor di et alJ 
l2007t> . 

(2) The amount of non-gravitational energy injection into 
proto-clusters and groups by AGN and starburst galaxies at 
intermediate-to-high redshifts z > 0.8 is poorly understood. 
Other observables are needed to constrain the physics and to 
answer this question which seems to be essential in under- 
standing the resulting gas profiles. Our simulations suggest 
that AGN-type feedback lowers the central pressure values as 
a hydrodynamic response of the gas distribution to the non- 
gravitational feedback of energy. This effect inhibits gas from 
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falling into the core regions which causes a flatter and more 
extended pressure profile and a noticeably reduced power of 
the SZ power spectrum at small angular scales for t > 2000. 

(3) For the SZ flux to be converged, an integration of the 
pressure profile out to 4/?200 is necessary; half of the SZ flux 
is contributed from regions outside 1?200- To compute a reli- 
able SZ power spectrum, it is essential to precisely charac- 
terize the state of the gas in these infall regions. In particu- 
lar, we find that: (i) the pressure support from kinetic energy 
strongly increases as a function of radius to reach on aver- 
age equipartition with the thermal energy at ~ 2/?200 in our 
AGN model with the exact dependence on cluster mass to be 
determined by future work; (ii) the effective adiabatic index 
T = dln/?/dlnp ~ 1.2 in the interior, but upturns towards 
r ~ 5/3 beyond the virial radius; (iii) the inclusion of cluster 
asphericity at large radii may also become important. 

Hence a successful semi-analytic model of the spherical 
cluster pressure, if that is indeed a viable goal, at the least 
needs careful calibration using numerical simulations which 
accurately treat all of the effects. The variance of the aver- 
age profiles also encodes important information that is man- 
ifested in the power spectrum. Our studies also show that 
simplified analytic models that employ hydrostatic gas mod- 
els with a constant F necessarily overpredict the SZ power on 
large scales by up to a factor of two and predict an inconsis- 
tent shape of the SZ power spectrum. The alternative that we 
explore in a subsequent paper is to use stacked scaled simu- 
lational clusters which are rotated to principal axes to provide 
the pressure form factors for the semi-analytic approach. 

The tSZ power spectrum of our 5 12 3 simulation agrees well 
with the average of our ten 256 3 simulations. A large number 
of simulations are needed to properly sample the high-mass 
end of the cluster mass function and hence accurately deal 
with sample (cosmic) variance. Alternatively, larger cosmo- 
logical volumes can compensate since they contain enough 
statistics on the large scale modes that are responsible in part 



for forming the highest-mass clusters which are also the rarest 
events. This, however, is quite challenging as we require the 
same (high-)resolution to accurately follow the physics in the 
cluster cores which is needed to obtain profiles that match 
current X-ray data. Our 256 3 simulations do not quite sample 
large enough scales to provide a fully converged kSZ power 
spectrum at low i since we miss the long-wavelength tail of 
the velocity power spectrum. We also have ignored the patchy 
re-ionization kSZ which could be a significant contr ibutor, up 
to 50% of the total kSZ (e.g.. lHiev et aLll2008ll2007l) . 

We have found the I < 2000 multipole range to be relatively 
insensitive to cooling and feedback, at least for the range con- 
strained by the X-ray data. We did find the higher multipole 
range (£ ~ 2000 - 10000) probed by the high-resolution ACT 
and SPT CMB telescopes is sensitive to the feedback prescrip- 
tion; hence the high-f SZ power spectrum can be used to con- 
strain the theory of intracluster gas, in particular for the highly 
uncertain redshifts > 0.8. In addition to the SZ power spec- 
trum probe, our simulations can be used to address the cosmo- 
logical significance of cluster counts as derived from the SZ 
effect. Counts provide complementary constraints on param- 
eters that help to break some degeneracies that are present in 
the power spectrum method. By employing inhomogeneous, 
localized and self -regulated feedback we are not only able to 
match recent X-ray reconstructions of cluster core regions, 
but also decrease the tension in <x 8 estimated from SZ power 
with <x 8 from other cosmological probes. However, only a 
detailed confrontation between simulations exploring the vast 
terrain of feedback options with the rapidly improving high 
resolution observations of cluster interiors can move the the- 
ory of cluster gas physics and its use for precision cosmology 
forward. 

We thank Norm Murray, Volker Springel, Hy Trac, 
Jerry Ostriker, Gil Holder, Niayesh Afshordi and Diasuke Na- 
gai for useful discussions. Research in Canada is supported by 
NSERC and CIFAR. Simulations were run on SCINET and 
CITAs Sunnyvale HPC clusters. 
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